import csv
import math
import random

# ========== 零依赖矩阵工具 ==========
def mat_shape(A):
    return (len(A), len(A[0]) if A else 0)

def mat_transpose(A):
    r, c = mat_shape(A)
    return [[A[i][j] for i in range(r)] for j in range(c)]

def mat_mul(A, B):
    ar, ac = mat_shape(A)
    br, bc = mat_shape(B)
    assert ac == br
    C = [[0.0] * bc for _ in range(ar)]
    for i in range(ar):
        for k in range(ac):
            aik = A[i][k]
            for j in range(bc):
                C[i][j] += aik * B[k][j]
    return C

def vec_sub(a, b):
    return [x - y for x, y in zip(a, b)]

def vec_dot(u, v):
    return sum(x * y for x, y in zip(u, v))

def vec_mean(vectors):
    n = len(vectors)
    dim = len(vectors[0])
    return [sum(vec[d] for vec in vectors) / n for d in range(dim)]

# ========== 特征值分解（幂法 + 放缩） ==========
def power_iteration(A, max_iter=1000, tol=1e-10):
    n = len(A)
    # 随机初始向量
    b = [random.random() for _ in range(n)]
    norm = math.sqrt(sum(x * x for x in b))
    b = [x / norm for x in b]

    for _ in range(max_iter):
        Ab = [sum(A[i][j] * b[j] for j in range(n)) for i in range(n)]
        lam = vec_dot(b, Ab)
        new_norm = math.sqrt(sum(x * x for x in Ab))
        if new_norm < tol:
            break
        b_new = [x / new_norm for x in Ab]
        if vec_dot(b, b_new) > 1 - tol:
            break
        b = b_new
    return lam, b

def deflate(A, lam, v):
    n = len(A)
    outer = [[lam * v[i] * v[j] for j in range(n)] for i in range(n)]
    return [[A[i][j] - outer[i][j] for j in range(n)] for i in range(n)]

def eig(A, k=None):
    if k is None:
        k = len(A)
    vals, vecs = [], []
    Ak = [row[:] for row in A]
    for _ in range(k):
        lam, v = power_iteration(Ak)
        vals.append(lam)
        vecs.append(v)
        Ak = deflate(Ak, lam, v)
    return vals, vecs

# ========== PCA 主流程 ==========
def center_matrix(X):
    col_means = vec_mean(X)
    return [vec_sub(row, col_means) for row in X], col_means

def cov_matrix(X):
    n = len(X)
    Xc, _ = center_matrix(X)
    Xt = mat_transpose(Xc)
    return mat_mul(Xt, Xc)  # Xt 形状 (d,n)，Xc 形状 (n,d) -> (d,d)
    # 上面返回的是 n * 协方差，因为没除 n；后续只关心方向，不影响

def pca(X, n_components=2):
    Xc, col_means = center_matrix(X)
    Sigma = cov_matrix(X)
    eigvals, eigvecs = eig(Sigma, k=n_components)
    # 按特征值从大到小排序
    order = sorted(range(len(eigvals)), key=lambda i: eigvals[i], reverse=True)
    eigvals = [eigvals[i] for i in order]
    eigvecs = [eigvecs[i] for i in order]
    # 投影
    W = mat_transpose(eigvecs)  # 形状 (d, n_components)
    Z = mat_mul(Xc, W)          # 形状 (n, n_components)
    return Z, eigvals, eigvecs, col_means

# ========== 加载数据 ==========
def load_csv(path, skip_cols=(0,)):
    """
    读取 CSV 并转成 float 矩阵。
    skip_cols: 需要忽略的列索引（例如 ID）
    """
    data = []
    with open(path, newline='') as f:
        reader = csv.reader(f)
        for row in reader:
            filtered = [v for i, v in enumerate(row) if i not in skip_cols]
            numeric = []
            for v in filtered:
                try:
                    numeric.append(float(v))
                except ValueError:
                    # 非数字字段（如职业、BMI 分类）简单哈希成数值
                    numeric.append(float(hash(v) % 100000) / 100000)
            data.append(numeric)
    return data

# ========== 运行 + 两位小数输出 ==========
if __name__ == "__main__":
    raw = load_csv("sleep.csv")   # 去掉 ID 列
    X = raw
    Z, lam, vec, mu = pca(X, n_components=2)
    print("睡眠质量数据集PCA：")

    # 特征值保留两位
    lam_round = [round(v, 2) for v in lam]
    print("前两个主成分方差:", lam_round)

    # 降维结果保留两位
    print("前 10 条样本降维结果:")
    for i, row in enumerate(Z[:10]):
        row_round = [round(x, 2) for x in row]
        print(i, row_round)